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Abstract 

We study accommodation of an excess charge carrier on long even-N polyynic oligomers 
CatH 2 due to displacements of the underlying carbon lattice and polarization of the surround- 
ing solvent in the context of carrier self-localization into a polaronic state. Spatial patterns 
of bond-length alternation, excess charge and spin densities are compared as derived with 
Hartree-Fock and two hybrid-DFT methods (BHandHLYP and B3LYP) in conjunction with 
the polarizable continuum model. Quite distinct resulting pictures of carrier accommodation 
are found when contributions from different interactions are analyzed. Solvation robustly acts 
to promote excess charge localization. 

Introduction 

A variety of one-dimensional (ID) semiconductor (SC) nanostructures such as ^-conjugated poly- 
mers (CPs), nanowires, nanotubes, and biological macromolecules attract a great deal of attention. 
They are interesting scientifically and can be exploited in various technological areas including 
(opto) electronics, energy harvesting and sensors. The nature and properties of excess charge 
carriers on these structures are fundamental for many processes and applications. For purposes 
of our discussion here, one distinguishes between nearly free band states of excess carriers and 
self-localized polaronic states. 1 ^ Self-localization and the formation of ID polarons may occur 
due to the (strong) interaction of the electronic subsystem with another subsystem such as dis- 
placements of the underlying atomic lattice, the mechanism extensively studied for various ID 
electron-phonon systems, and particularly for CPs. 5 " 7 Another, and less explored, implementa- 
tion of the polaronic effect can take place for ID SCs immersed in 3D polar media, the situation 
common for applications involving fundamental redox processes in polar solvents. In this case of 
what could be called charge carrier solvation, the long-range Coulomb interaction is expected- - — 
to lead to the formation of a localized electronic state on a SC structure surrounded by a self- 
consistent pattern of the sluggish (orientational) polarization of the solvent. Among consequences 
of self-localization, lattice-deformation and solvation-induced polarons can feature new signatures 
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in the optical absorption due to local intragap electronic levels- 2 ^""-^ 2 -^ an( j substantially reduced 
mobilities ,2JLl£.dl 

In this paper we use three different theoretical levels of ab initio computations to compare 
pictures of accommodation of an excess electron on a prototypical ID SC system: long even-Af 
polyynic chains C^H2, emerging due to both bond-length adjustments of the carbon lattice and 
polarization of the surrounding solvent. The system we study is nearly charge-conjugation sym- 
metric and similar results are obtained 11 upon an addition of an extra hole. In even-/V oligomers 
the arising patterns of bond-length modulations are non-topological (polaron-like) being thus rep- 
resentative of a broad range of CPs. A related comparative study of odd-N polyynic oligomers 
has been recently published; 16 those feature topological kink-like patterns that can be found in 
polymers with the degenerate ground state, trans-polyacetylene being another example.- - - 

In applying first-principles calculations one expects to be able to elucidate effects due to reor- 
ganization of valence electrons and many-electron interactions, which can be crucially important in 
accommodation of excess charge carriers. Numerous studies at various theoretical levels have been 
published on the structure of the ground state and excitations in CPs. It is therefore noteworthy that, 
despite a relatively long history of electron-lattice polarons in CPs (see, e.g., multiple references in 
Refs.— ^), some questions have been raised recently regarding applicability of different ab initio 
frameworks to describe the polaron formation. That pertains to the reported failure of the local- 
density-approximation and generalized-gradient-approximation DFT (density-functional theory) 
schemes to detect self-localized charge density distributions in various charged oligomers, while 
Hartee-Fock (HF), parameterized semi-empirical and possibly hybrid-functional DFT 18 methods 
have been reported to lead to charge localization in the middle of oligomers Quite dif- 

ferent level-of-theory-dependent magnitudes of the effective electron-phonon coupling have also 
been found for the dimerized ground states of polymeric systems.-^ 2 -" 2 ^ These observations thus 
bear on an important general issue of a choice of appropriate ab initio methods 18,24 to faithfully 
describe properties of ID semiconductors. In our computations for this paper, we compare results 
derived with the GAUSSIAN 03 implemented 25 HF and two popular hybrid-functional DFT meth- 
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ods (BHandHLYP and B3LYP) in vacuum and solvent environments, the latter represented within 
the framework of the polarizable continuum model.— >^ We will be interested in the carbon-carbon 
bond-length alternation (BLA) patterns and spatial distributions of excess charge and spin densities 
resulting from accommodation of an extra electron. 

Our choice— dUS. of linear polyynic carbon chains as a model system is to a large extent related 
to their structural simplicity that allows us to keep computations on relatively long oligomers (up 
to N = 100 in this paper) still practicable. In doing so, we attempt to address the intrinsic response 
of the system as determined by inherent interactions rather than by the end effects due to the hydro- 
gen termination. In terms of the polaron formation, one would like to see if the polaron-like spatial 
structure is indeed self-maintained, that is, if its spatial extent is not affected by the length of the 
oligomer. We have explicitly demonstrated 16 in the case of kink- soli tons that their resulting spatial 
extent strongly depends on the computational method used. Our results in this paper also display a 
very strong method-dependence as will be discussed later. While being a model system, one should 
note that polyynic chains continue to be the subject of much attention in their own hght^^^ r — 
A qualitative difference in the electronic structures of polyyne and trans-polyacetylene is that the 
former features an extra spatial degree of freedom: ^-electron molecular orbitals can have two 
independent orientations perpendicular to the polymer axis. Within the framework of appropri- 
ately modified Su-Schrieffer-Heeger (SSH) model of conjugated polymers,-- it was show n 34 ' 
that, thanks to the extra degeneracy of ^-electron levels, polyynes would possess a particularly 
rich family of electron-lattice self-localized excitations. Also, spin-charge relationships for such 
excitations are quite different from trans-polyacetylene's. 

In order to distinguish between various effects (due to electron-electron interactions, solvent 
polarization and carbon lattice adjustments) on resulting charge and spin densities, we study both 
systems with prescribed carbon atom positions (referred to as "rigid geometry" cases) as well 
as fully optimized systems, in which carbon atom positions minimize the total energy. Rigid 
geometries feature uniform BLA patterns chosen based on the results of optimization in the ground 
state of neutral oligomers Ca^Hj. The data on neutral oligomers is used as benchmarks against 
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which we determine changes in BLA patterns and excess charge distributions in negatively charged 
oligomers C^H 2 . 

Computations and data processing 

All ab initio computations in this study were performed using the GAUSSIAN 03 suite of pro- 
grams. 25 We have employed pure Hartree-Fock (HF) along with the well-known hybrid-DFT 
functionals BHandHLYP and B3LYP in order to compare results of the different levels of the- 
ory. Hybrid exchange-correlation density functionals 18 include both local and non-local effects 
and are commonly used in studies of electronic properties of CPs. The reader is reminded that 
functionals used here differ by the amounts of HF exchange included: 50% for BHandHLYP and 
20% for B3LYP. A rich 6-31 l++G(d,p) all electron basis set has been employed throughout. To 
investigate the effects of solvation in a polar medium, GAUSSIAN 03 25 offers its implementation 
of the polarizable continuum model (PCM) described in original publications.-^"— For all "in sol- 
vent" PCM calculations, water has been chosen as a polar solvent with its default parameters in 
GAUSSIAN 03. In working with long oligomers we met with certain method-dependent limita- 
tions on obtaining satisfactory levels of convergence. In order to produce reliable results within the 
HF method, restricted (R) and restricted open (RO) shell wave functions were used. For the hybrid 
DFT methods BHandHLYP and B3LYP fully unrestricted (U) wave functions were utilized. Note 
that, within each computational method, we consistently used the same form of the wave function. 

The subject of our attention in this paper are long QyH2 chains with even number of carbons 
N (N up to 100). It is well established 22 ^^^^ that long even-N neutral oligomers have their 
ground state with the polyynic structure, that is, they feature an alternating pattern of triple C = C 
and single C — C bonds and a gap in the electronic spectrum. We note that wherever possible, we 
verified that results of our "in vacuum" calculations compare well against previously publishe d 23i3 ° 
ab initio data on even-Af neutral and charged QyH2 systems both in terms of energetics and opti- 
mized bond lengths. To our knowledge, oligomer lengths explored in our studies (Refs.— ^ and 
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this paper) are appreciably longer than used previously in comparable calculations for charged 
systems. 

A very convenient and well-known- way to characterize the geometry of dimerized poly- 
meric structures is via bond-length alternation (BLA) patterns 

8 n = (-l) n (l n -l n - l ), (1) 

where l n is the length of the nth carbon-carbon bond defined, e. g., as on the right of the nth 
carbon atom (we do not discuss the end bonds to hydrogens). In the infinite dimerized neutral 
structure, the dimerization pattern ?? would be uniform, that is, independent of the spatial bond 
position n, we denote the magnitude of this equilibrium ground-state pattern as 8q (this is not to 
be confused with a specific site value of 5 n= o). The double degeneracy of the ground state of the 
infinite polymer is in that the uniform pattern ?? can assume either value of +6q or — 6q. 5 ~ 7 For our 
purposes, we retrieve the magnitude 6q from the middle segment of long even-N neutral oligomers 



featuring a well-established spatially-independent behavior, see illustration in [figure][l][]l , and 



display that ground-state BLA pattern as positive. As mentioned earlier, the magnitude Sq is found 



to be strongly method-dependent; this dependence is reflected in results shown in [figure] [3] []3 
and in published 16 ' 23 tabulated numerical data. 

While exhibiting a spatially-independent behavior in the middle of long neutral oligomers, 



[figure] [1][] 1 'a) also clearly shows that the end (oligomer termination) effects on the BLA pattern 
are quite substantial. In numerous computations we observed that the "shape" of the end-specific 
BLA structures stabilizes and becomes easily discernible in long enough oligomers (of course, 
within a specific computational method). Practically the same persistent shape has been found 
in end BLAs of neutral and charged, even- and odd-N oligomers This fact allows to "elimi- 
nate" the end effects by comparison and an appropriate "subtraction" procedure, as illustrated in 



[figure] [1][]1 b). Such a procedure of removing the end effects from the finite oligomer data has 



been used to produced the BLA patterns shown in [figure] [3][]3 focussing thus on the changes 
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Figure 1: Illustration of end effects in BLA patterns and their elimination: Panel (a) compares 
BLA patterns for the neutral (0) and charged (-1) chains with end effects. Panel (b) indicates the 
BLA pattern after the subtraction of common end BLA structures. Vacuum results are shown with 
dashed lines and solid lines indicate the resulting BLA pattern upon solvation. Results depicted 
here were obtained on the C100H2 chain using the B3LYP method. 



induced by the accommodation of an excess charge carrier. 

Very essential for our discussion are the resulting spatial charge and spin distributions over the 
oligomers, which are naturally related to atomic charges and spins in outputs of ab initio calcula- 
tions. Importantly, these quantities are calculated from many-electron wave functions and thereby 
reflect responses due to all electrons in the system. Two procedures, after Mulliken and Lowdin, 
are widely used for charge population analysis 26 and only Mulliken procedure is available for 
spins in GAUSSIAN 03. It is well known 26 that calculations of atomic-centric quantities depend 
on the basis set transformations and sometimes lead to artifacts and spurious results. We have 
demonstrated previously 11 that certain cancelation effects take place and Mulliken and Lowdin 
procedures lead to more similar results when one is interested in the spatial distribution of the 
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Figure 2: Illustrations of the spatial behavior of raw atomic densities scaled with respect to their 
maximum values. Both examples are derived from oligomers with uniform rigid geometries, (a) 
Excess charge density obtained for C^qH^, m vacuum with the HF method, (b) Spin density ob- 
tained for C100H2 in solvent with the BHandHLYP method. 



excess charge derived as the difference of charge densities on charged and neutral chains. An 
example of Lowdin-derived excess atomic charge densities q n from our current computations is 



shown in [figure] [2] []2 ^a). [figure] [2] []2 ^b) illustrates an instance of Mulliken-derived atomic spin 
densities p;,. When used for absolute data values, we keep both densities normalized to unity: 



£<7n = Y,Pn 



1. 



In [figure] [2] [] 2, the data is shown in a scaled form to emphasize its spatial behavior. It is evident 
from this figure that atomic-centric quantities generally exhibit very strong variations from carbon 
to the next carbon. It is not clear to what extent these oscillations may result from computational 
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artifacts perhaps also related to the presence of bond- oriented density waves- - - in our system. In 
order to eliminate (possibly spurious) high spatial-frequency effects and to focus on longer spatial- 
range structures, some averaging/smoothing of the atomic data can be applied. Needless to say 
that results can thus depend on details of the smoothing procedure. In our earlier work,— we have 
shown that a simple cell averaging of excess charge data can already lead to quite smooth results, 
moreover it has made Mulliken- and Lowdin-derived results nearly identical. The unit cells here 
are defined as consisting of pairs of neighboring carbon atoms starting from the oligomer ends (as 
is common, the chain end cells also include hydrogen atomic charges/spins). We apply the same 
summation over the neighboring carbons to generate cell-centric data in this paper. To smooth 
data even further, we also apply a simple-minded nearest-neighbor-cell averaging where a cell- 
specific quantity Xj is transformed into (xj) = 0.5xi + 0.25 (.*,■_ i +x,- + i). So generated smoothed 



data on cell-centric charge and spin p- densities are used in [figure] [3] []3 While the number 
of cells (index /) is of course twice as small as the number of carbons, the procedure preserves the 
normalization: 

X>=X>' = 1 - 

i i 

While reiterating the caveat on the procedure dependence (especially on some spin data), we men- 
tion here that the described above averaging produce results with smaller-amplitude variations 
than some other approaches we tried, which is a desirable feature in the absence of other criteria 
to ascertain the physical validity of the atomic-centric data from computations. 

As discussed in the Introduction, in addition to the fully relaxed geometries of the underlying 
atomic lattice that minimize the total system energy, we also study rigid geometries (RG) with pre- 
scribed bond lengths. Those have been generated separately for each of the computational methods 
used in respective environments (vacuum and solvent). We would first find a fully optimized ge- 
ometry of the N =100 neutral oligomer and retrieve the resulting lengths of the C = C and C — C 
bonds in its central part (as well as the length of the terminating C — H bonds). We then assign 
those bond lengths alternatingly to all neighboring carbon pairs throughout the oligomer. The re- 
sulting RG structures thus exhibit entirely spatially uniform carbon-carbon dimerization patterns 



within the confines of terminating hydrogens. We note that there is very little difference, as clearly 
seen in | [figure] [!][]![ in the dimerization amplitude do between neutral chains in vacuum and in 
solvent. 



Results and discussion 

In this Section we discuss our computational results processed as described above. Our discussion 



will revolve around [figure] [3] [] 3 displaying spatial distributions of BLA patterns, excess charge 
and spin densities for N = 60 (ROHF computations) and N = 100 (UBHandHLYP and UB3LYP) 
polyynic oligomers C^H^ with an extra electron. 

It is instructive to first recall how self-trapping of an excess charge carrier (an electron or a 
hole) is described in a single (quasi)particle picture within a standard ID continuum adiabatic 
framework. 3 The self-trapped state corresponds to the localized ground-state wave function y(x) 
in the solution of the non-linear Schrodinger equation: 



h l d 2 y{x) 
2m dx 2 



- Jdx'V (x-x')\ y/-(x') 1 2 yf(x) = E yf(x) , (2) 



where x is the coordinate along the structure axis, m the intrinsic effective mass of the carrier, and 
V(x) the effective self-interaction mediated by another subsystem. In the case of a short-range 
electron-phonon mediation, the self-interaction can be taken local: V(x) = gd(x) leading then 
to the widely known exact result— 5 ^ 3 Y(x) 1/ cosh(gmx/2h 2 ) for a continuum ID polaron. In 
the case of the long-range polarization interaction, the effective V(x) behaves as l/\x\ at large 
distances, while the short-range behavior depends on specific system details such as, e.g., the 
actual transverse charge density distribution and the geometry of the dielectric screening.-^ In 
the absence of the self-interaction, ?? would just yield spatially delocalized plane- wave solutions. 
Among features of this single-particle description is that the same probability density |i//(jc)| 2 de- 
termines both spatial distributions of charge and spin of the charge carrier. The other one is that 
it involves the adiabatic framework, that is, the mediating subsystem is assumed very slow allow- 
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Figure 3: Processed results for charged C^H 2 oligomers as arranged in columns by the compu- 
tational method: (a) ROHF, (b) UBHandHLYP, and (c) UB3LYP. Row 1 displays optimized BLA 
patterns with end effects eliminated for oligomers in vacuum (dashed lines) and in solvent (solid). 
The notation is different in rows 2-5, where solid lines show results for fully optimized carbon 
lattice geometries, while dashed lines (when distinguishable) for rigid geometries with uniform 
dimerization patterns (see the text). These rows display cell-centric excess charge and spin densi- 
ties for chains in vacuum (rows 2 and 3) and in the solvent environment (rows 4 and 5). 
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ing to determine electronic states for its given static configuration, which in turn self-consistently 
responds to the probability | y(x)\ 2 ("strong-coupling polarons"- - -). It should be emphasized that 
in realizations of the computational methods used, both displacements of the underlying carbon 
lattice and solvent polarization constitute such slow subsystems. 

The single-particle description has the drawback of not explicitly including valence band elec- 
trons, whose reorganization may be important in the accommodation of the excess carrier. Such a 
reorganization can take place even in electron-lattice models without direct electron-electron inter- 
actions that have been used to describe non-linear localized excitations in CPs. 5,6 A nice analysis^ 
was in fact given of how polarons of a two-band Peierls dielectric model of CPs evolve into single- 
particle Holstein polarons in the limit of the "frozen valence band" approximation. Polaronic and 
solitonic excitations in these models are also derived within self-consistent adiabatic treatment of 
the underlying lattice. In ab initio computations, reorganization of valence band electrons oc- 
cur as driven by both the interaction with the bond-length modulations and by electron-electron 
Coulomb interactions. With many electrons determining the resulting spatial densities, excess 
charge and spin distributions should not be necessarily the same. One of the famous instances 
of this phenomenon is inverse charge-spin relationships for kink-solitons in trans-polyacetylene- 5 -^ 
(these relationships are different for kinks in polyynes— >2 5 -). We will also mention in this regard 
separation of charge and spin degrees of freedoms well known in ID electronic liquids. 45,46 Still 
one could expect that, for non-topological localized polaronic excitations, excess charge and spin 
density distributions would closely follow each other. 

Turning now to results in [figure] [3][]3[ we note upfront that, similarly to our observations on 
kink-solitons, 16 spatial extents of the structures developed as a result of electron accommodation 
are strongly computational method-dependent. The magnitudes 8q of the uniform dimerization are 
consistent with the notion of the increased effective electron-lattice interactions with the increase 
of the amount of the HF exchange in the method functionals. We would be more interested now in 
qualitative aspects of the results as related to the phenomenon of self-localization. 

The main reason for our studying rigid geometry (RG) cases has been to try to distinguish ef- 
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fects due to different interactions. As RGs feature spatially uniform BLA backgrounds, one might 
expect that in vacuum the excess electron would be delocalized over the oligomer showing the 
particle-in-the-box type behavior for charge and spin densities. On the other hand, if the lattice 
is allowed to relax, that could cause, as per a conventional adiabatic polaron picture, spatial lo- 
calization of excess charge and spin. Vacuum results for the densities are shown in Rows 2 and 



3 of [figure] [3] [] 3 One indeed observes that RG data displays a delocalized behavior for both 
charge and spin when calculated with hybrid-DFT methods. Not so however when calculated in 
HF. Hartree-Fock data (panels (a2) and (a3)) display localization of excess charge and spin already 
in RG, and the degree of localization is only very weakly affected by a quite substantial dip in 
the optimized BLA pattern (panel (al)) when the lattice is allowed to relax. We conclude that 
self-localization observed in these results is not caused by the lattice deformation but rather by 
electron-electron interactions as treated in HF. The deformation of the BLA pattern then results 
as a response to the developed localized density. We easily confirmed this kind of localization 
due to electron-electron interactions alone in simplified calculations once these interactions are 
treated in the HF approximation. One cannot help wondering if this qualitative picture may be 
just an artifact of the HF treatment. It is not clear if there are slow (low-frequency) electronic 
excitations in this system that exhibits a well-known very substantial HF energy gap.— ^ An al- 
ternative to that picture would be that as a result of Coulomb interactions with the fast electronic 
subsystem, the excess charge carrier remains delocalized but with its properties renormalized (see, 
e.g., calculations of corrections to the DFT quasi-particle energies within the framework of the 
GW-approximation-^l^). Such a renormalized quasi-particle could then interact with the lat- 
tice displacements to form an electron-lattice polaron. Perhaps a properly realized time-dependent 
computational scheme could clarify the question of the purely electronic response. 

While yielding a delocalized behavior in RG, hybrid-DFT computations exhibit other issues in 
the context of the polaron formation due to carbon lattice displacements in vacuum. Particularly, 
as has already been noticed,— we do not observe electron self-localization in results obtained 
with UB3LYP - the slight changes in densities due to relaxation (panels (c2) and (c3)) and a 



13 



flat optimized BLA pattern in the middle of oligomer (panel (cl)) indicate that they are finite 
oligomer-size effects - at least within the range of oligomer lengths studied. With UBHandHLYP 
data, visually one could probably identify the excess charge density (panel (b2)) in the optimized 
system as localized. There is however no counterpart to this behavior in the spin density (panel 
(b3)), which rather looks like a spin-density wave (SDW). Also quite unusual is the shape of the 
dip in the optimized BLA pattern in panel (bl), quite substantial in magnitude but the spatial extent 
of which is still likely limited by the oligomer ends rather than being self-consistent. 

Concluding the survey of the vacuum data in rows 2 and 3, none of the methods we used seems 
to have produced results that can be qualified with certainty as self-localization of the excess elec- 
tron charge and spin densities caused by its interaction with the displacements of the underlying 
carbon lattice. 

Quite substantial effects take place due to the polarization of the surrounding solvent by the 



excess electron. They are easily discernible in [figure] [3] [] 3 One immediately notices in row 1 
that optimized BLA patterns in the solvent acquire a (much) more spatially localized character. 
Results for the charge and spin densities in the solvent environment are shown in rows 4 and 5. 

Using B3LYP computations, we have previously demonstrated 11 that the solvation effect alone 
can lead to self-localization of excess charge densities. RG data in panels (c4) and (c5) clearly show 
that this is the case for both charge and spin densities, even though a complete self-consistency is 
not completely achieved yet at this oligomer length. The behavior of spin and charge distributions 
here closely follow each other. As HF computations yield localized patterns already in RG, one 
cannot study the solvation effect here separately. It is however evident from panels in column 
(a) that solvation results in increased spatial localization of excess charge and spin (all HF data 
exhibit a noticeably larger degree of spin localization than that of charge). The BHandHLYP data 
in column (b) also shows a significant degree of localization of charge density due to solvation 
(b4), appreciably more than what is caused by the lattice deformation in vacuum (b2). Once again, 
though, the spin response (b5) in these computations for RG is quite surprising: while, perhaps, 
exhibiting a localized pattern, it is accompanied by what appears to be a large amplitude and spatial 
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scale SDW-like response. One, in fact, would not easily discern a localized pattern in the raw 



atomic data shown in |[figure][2][]2| ^b), and different averaging procedures could lead to different 
conclusions. The spin density distribution in this case however does acquire a more delineated 
localized shape, even if on the background of a smaller SDW-like pattern, when the carbon lattice 
is allowed to relax. The lattice relaxation does not appear to affect the already well-localized excess 
charge density in panel (b4). Very small, and apparently in different directions for charge and spin, 
is the lattice relaxation effect in the HF data ((a4) and (a5)). The only computational method 
that clearly shows a synergistic amplifying localization effect of solvation and lattice relaxation 
on both spin and charge densities is B3LYP ((c4) and (c5)). One can however say that when 
solvation and lattice relaxation act together, they do cause the localization of excess charge and 
spin densities in BHandHLYP computations as well (as compared to the delocalized behavior 
in the vacuum RG environment). Overall, it is evident that, independently of the computational 
method, solvation always acts as to boost the degree of localization of excess charge density - the 
result quite understandable on simple physical grounds. It is hard to say if the peculiarities of the 
spin responses (and perhaps related BLA patterns) found in BHandHLYP results are caused by 
some computational artifacts or may indeed have a physical origin in different behaviors of charge 
and spin degrees of freedom due to many-electron interactions. It would be interesting to see if 
similar results are obtained in long C4„ + 2H 2 oligomers that have a different electronic symmetry— 
than Cz^Hj systems studied in this paper. 

A truly self-consistent localized solution should be (practically) oligomer-length independent, 
which could be straightforwardly verified by comparing all spatial patterns (BLA, charge and spin 
density distributions) on oligomers of increasing lengths. That is, for instance, the case when we 
compare relevant results obtained in the HF computations for N = 40 and N = 60 oligomers. On 
the other hand, vacuum BLA patterns obtained in hybrid-DFT computations appear to continue 
to evolve even at our current limit of N = 100 oligomers. The difference between results derived 
for oligomers of different lengths is substantially much smaller for computations in the solvent 



environment. This is in agreement with the observation evident from [figure] [3] [] 3 that solvation 
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generally leads to more spatially localized patterns, which are therefore much weaker affected by 
chain ends in our longest oligomers. 



One can also look at the self-consistency issue from the viewpoint of energetics. |[table][l][]T 



lists some energy parameters derived from the computational outputs for oligomers of different 
lengths. A well-defined and physically important quantity characterizing accommodation of an 
excess electron is the electron affinity Ea, corresponding to the difference of total system energies 



for neutral C^H2 and negatively charged C^H^ chains. | [table] [ 1 ] [] 1 1 features electron affinities 
El in vacuum and £| in solvent environments, all obtained for fully optimized systems. Needless 
to note that numerical values in the table are strongly method-dependent. One can easily see 
that, within each method, the values of £| change only little when going to longer oligomers, an 
indication of the convergence to self-consistent solutions. Vacuum values E^ change by larger 
quantities, which is contributed to by two factors: vacuum structures are farther away from self- 
consistency, and the unscreened Coulomb field of the excess charge is long-range (long-range 
electron polarization effects may be reflected in the illustration of [figure] [2] []2 ^a), see also Ref. 16 ). 
The solvent environment efficiently screens such long-range Coulomb effects.— 
Table 1: Electron affinities and polarization energies in eV 



method 


system 


EX El U° U 


ROHF 


C40H2 
C60H2 


0.88 2.10 0.38 1.84 
0.92 2.11 0.44 1.87 


UBHandHLYP 


C60H2 
C100H2 


3.00 3.82 0.34 1.39 
3.11 3.85 0.40 1.40 


UB3LYP 


C80H2 
C100H2 


3.92 4.41 0.30 0.94 
4.00 4.42 0.32 0.92 



|[table][l][]l| also shows electrostatic energies U° and U stored in the solvent polarization 
for, respectively, neutral and charged oligomers, as retrieved from GAUSSIAN 03 outputs. We 
previously 11 analyzed the evolution of these quantities with oligomer lengths and found that U° 
grows linearly with the length while U first decreases with the length, achieves its minimum 
when the length becomes larger than the self-consistent size of the self-localized charge carrier, 
and then continues to increase just as U° does. By these features, one can also judge from the table 
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data that self-consistency has been practically achieved in ROHF and UBHandHLYP computations 
but longer oligomers would need to be employed to observe explicit self-consistency in UB3LYP 
computations. 

Looking at the trends displayed in results obtained with two hybrid-DFT methods, we cannot 
exclude that there may be a range of "customized" hybrid functionals with the "right" mixture of 
HF and DFT correlations that would yield results featuring self-localization behavior for an excess 
charge carrier consistent with the single-particle picture of ??. One could then observe that in the 
absence of the interaction with a slow subsystem (uniform RG in vacuum) the resulting charge and 
spin densities are delocalized, while self-localization takes place due to the lattice deformation and 
due to solvation as separate effects as well as synergistically. In addition, the resulting localized 
charge and spin densities would largely follow each other even if with small variations due to many- 
electron interactions. Another conjecture also suggested by our data is that the failure to detect the 
formation of electron-lattice polarons in pure DFT and B3LYP computations may be related to the 
unscreened long-range Coulomb field of the excess charge in vacuum. Indeed, we clearly detected 
the formation of lattice kink-solitons in such computations, 16 and electron-lattice polarons can 
be thought of as bound states of kink-antikink pairs. 5 - 6 The screening of long-range Coulomb 
forces by fast polarization of the surrounding medium could energetically stabilize the polaron 
formation. Fast polarization here would act only as to reduce the magnitude of Coulomb forces. 
This screening effect is different from nearly static reorganization of a slow solvent considered in 
this paper and demonstrated to be a robust promoter of excess charge localization. Exploring the 
above conjectures could be an interesting avenue for further developments. 
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